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ABSTRACT 

We study the r-modes propagating in steadily mass accreting, nuclear burning, and 
geometrically thin envelopes on the surface of rotating neutron stars. For the modal analysis, 
we construct the envelope models which are fully radiaitive or have a convective region. We 
simply call the former radiative models and the latter convective models in this paper. As the 
angular rotation frequency fl is increased, the oscillation frequency u> of the r-modes in the 
thin envelopes deviates appreciably from the asymptotic frequency ui = 2m£l/l'{V + 1) defined 
in the limit of f2 — > 0, where u> is the frequency observed in the corotating frame of the star, 
and m and V are the indices of the spherical harmonic function Y™ representing the angular 
dependence of the modes. We find that the amplitudes of the fundamental r-modes with no 
radial nodes of the eigenfunctions are strongly confined to the equatorial region, and to becomes 
only weakly dependent on f2, gathering in a frequency range of lo/2tt <10Hz, at rapid rotation 
rates. We also find that the fundamental r-modes in the convective models are destabilized 
by strong nuclear burning in the convective region. Because of excessive heating by nuclear 
buring, the corotating-frame oscillation frequency u of the r-modes in the convective models 
becomes larger, and hence the inertial-frame oscillation frequency |er| becomes smaller, than 
those of the corresopnding r-modes in the radiative models, where a — u> — mfl is negative for 
the r-modes of positive m. We find that the relative frequency change / = — (a C onv — <Jrad)i 1 o ra d 
is always positive and becomes less than ^0.01 for the fundamental r-modes of V > \m\ + 1 at 
Wrad\/2Ti ^300Hz for m = 1 or at |<7 rac ;|/27r ~600Hz for m = 2, and that we need to consider 
the r-modes of I' much larger than \m\ for values of / as small as ^0.001, where a conv and 
(J ra d denote the oscillation frequencies for the convective and the radiative envelope models, 
respectively. 
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1. Introduction 

X-ray flux fluctuations observed in type I X-ray bursts in LMXB's (see, e.g., Strohmayer et al 
1996, 1997 for so called burst oscillations) are quasi-periodic oscillations whose frequencies are found in 
a frequency range around ^300Hz or around ~600Hz. During a burst, from the onset to the tail, the 
oscillation frequency in a given source usually increases with time by a small amount and gets saturated 
to an asymptotic frequency in the tail. The observed frequency drifts are generally less than ~1% of the 
mean frequency (e.g., Cumming & Bildsten 2000). However, the time evolution of burst oscillations are not 
always that simple. In fact, some X-ray burst sources have been reported to show a temporal decrease of 
the oscillation frequency in the tail of bursts (e.g., Strohmayer 1999, Miller 2000, Muno et al. 2000), or to 
show simultaneously two oscillation frequencies separated by a few Hz in a burst (Miller 2000, Galloway 
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ct al 2001, Muno et al 2002). Burst oscillations could be very important for neutron star physics. On the 
assumption that the oscillation frequencies are directly tied to the spin frequency of the star, the stability 
of the asymptotic frequency found in the tail of bursts separated by a few years has been used to determine 
some of the physical parameters of binary systems (e.g., Strohmaycr et al 1998, Gilles et al 2002). For 
recent observational and theoretical developements concerning burst oscillations as well as thermonuclear 
bursts on neutron stars, see a review paper by Strohmayer & Bildstcn (2003). 

As a model for burst oscillations, Strohmayer et al (1997) proposed that the oscillations are produced 
by inhomogeneous hot burning spots, the burning fronts of which are expanding laterally on the surface 
of a rotating neutron star, and that the frequency drifts are associated with the vertical expansion and 
contraction of the burning spots conserving angular momentum. The model was later refined by Cumming 
& Bildsten (2000), and Cumming et al (2001) concluded that the possible amount of relative frequency 
changes expected in the model was by a factor of 2 to 3 smaller than the largest observed values (see also 
Muno, Ozel, & Chakrabarty 2002). This model also has a difficulty in explaining the persistence of the 
oscillations to the tail of a burst, in which no inhomogeneity to produce the oscillations will remain because 
all the surface area has been reached by the burning fronts. Spitkovsky, Levin, & Ushomirsky (2002), on the 
other hand, discussed the global hydrodynamical flow in the ocean of an accreting neutron star, by taking 
into account the rapid rotation and the lift-up of the burning ocean, and conjectured that Jupiter-type 
vortices created in strong zonal fluid currents driven in a nuclear burst are responsible for burst oscillations 
in the tail. 

It was Heyl (2001) who suggested that the burst oscillations could be produced by the r-modes of low 
azimuthal wavenumber m = 1 or m = 2 traveling in the surface ocean of the neutron star. The r-modes are 
retrograde toroidal modes and are called Rossby waves in geophysical context (see, e.g., Greenspan 1968, 
Pedlosky 1987). The r-modes in mass-accreting envelopes on the surface of a slowly rotating neuton star 
have been discussed by Strohmayer & Lee (1996), who found that the fundamental r-modes with no radial 
nodes of the eigenfunctions are all pulsationally stable except for the cases of mass-acceretion rates M as 
low as M/MEdd ~ 0.001, where MEdd denotes the Eddington rate. In this paper, we extend the previous 
analysis to the cases of rapidly rotating neutron stars, considering that the spin frequencies fi of many 
neutron stars in LMXB's are now believed to be in a rather narrow range around fl/2-K <~ 300Hz. The 
main points behind the r-mode model for burst oscillations may be described as follows. Let us assume 
that strong nuclear burning in the convective phase of a nuclear burst excites surface r-modes whose time 
and azimuthal dependence is given by exp(i[u> conv t + m<p]), where oj denotes the oscillation frequency in 
the corotating frame and the integer m is the azimuthal wave number. For the r-modes, ui is proportional 
to and is much smaller than the spin frequency f2 of the star, and u) is positive for positive m in the 
convention we use in this paper. Observers measure the inertial frame frequency o~ conv = Lu conv — mCt, 
the magnitude of which is |<r cont ,| — inVt — ^conv- 

At the end of the burst, the ocean cools and becomes 
radiative, and the rotating frame oscillation frequency decreases to ui ra d and the magnitude of the inertial 
frame frequency increases to |<7 ra£ 2| = rreO — u> ra d- This model then predicts a fractional increase in the 
measured inertial frame frequency from the start of the burst to the asymptotic end, which may be given by 
/ = —(cr conv — a r ad) I o~rad- To estimate possible magnitudes of the quantity /, we calculate in this paper the 
r-modes for radiative envelopes and convective envelopes, assuming that the envelopes are in steady state. 
§2 and 3 briefly describe method of calculation employed in this paper, and §4 is for numerical results, and 
§5 and 6 are for discussion and conclusion. 
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2. Envelope Calculation 

Method of calculation for steadily mass-accreting, nuclear-burning, thin surface envelopes is almost the 
same as that given in Strohmayer & Lee (1996). No effects of general relativity are included in this paper. 
If we use the independent variable defined as 

z = M-M ri (1) 

where M and M r are the total mass of the star and the mass contained in the sphere of radius r, the basic 
ordinary differential equations for mass-accreting and nuclear-burning thin envelopes in steady state are 

* 1 (2) 



dz Aitr 2 p' 

dp _ GM r 
dz Airr 4 ' 



(3) 



dL r • f dU p <ip\ 

= -e N + e„ + M — - 4-f , (4) 



dz " \dz p 2 dz 

and the energy transfer equation in convection zones is 

dT GM r T ( <91nT\ GM r T 



dz 4irr 4 p \d\np J ad 4irr 4 p^ ad ' ^ 
and that in radiative zones is 

dT 3KL r . . 

= 647r 2 acr 4 T 3 ' ( ' 

where p is the mass density, p is the pressure, T is the temperature, U is the internal energy per gram, k is 
the opacity, 6n is the nuclear energy generation rate per gram, e v is the neutrino energy loss rate per gram, 
L r is the luminosity at r, M = —A-Kr 2 pv r is the mass accretion rate, and G, a, and c are the gravitational 
constant, the radiation constant, and the velocity of light, respectively. The structure of the chemical 
composition, X i: in steady state is governed by 

where rrij is the mass of the particle j, Ti is the rate of the reaction i, and riji is the number of particles j 
created (nji > 0) or destroyed (riji < 0) in a reaction i. The equation of state p = p(p,T,Xj) is that of a 
perfect gas with radiation, and is given by p = pi + p e + p r , where pj — niksT with nj and ks being the 
number density of ions and the Boltzmann constant, p e is the electron gas pressure, and p r = aT 4 /3 is the 
radiation pressure. The electron gas can be arbitrarily degenerate and relativistic. For hydrogen burning 
energy generation rates, we use the formulae given in Reeves (1965). The reaction rates of the CNO cycle 
in high temperature environment are limited by the /3-decays of 14 and 15 0, leading to the saturation of 
the rates to t^% x = 5.9 x lQ lb Z CN o (e.g., Wallace & Woosley 1981), where Z C no = 12 C + 14 N + 16 O. 
For helium burning we include the 3a and 12 C(a,-f) 16 reactions, the rates of which are given in Fowler, 
Caughlan, & Zimmerman (1975). The rate of the carbon burning 12 C + 12 C reaction is from Patterson, 
Winkler, & Zaidins (1969), and Arnett & Truran (1969). We apply both weak and strong electron screening 
corrections to the nuclear energy generation rates following Salpeter & Van Horn (1969). For the energy 
loss rate due to neutrino emission we employ formulae given by Beaudet, Petrosian, & Salpeter (1967). For 
opacity, we use a subroutine written by Paczyhski (1970) who adopted Cox & Stewart (1970) radiative 
opacity and the conductive opacity by Hubbard & Lampe (1969) and Canuto (1970). 
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For given M, R, and M, we integrate the set of differential equations from the surface z = to the 
bottom z — z of the envelope. At the surface of the envelope, where we assume r = R and M r — M, the 
surface tempetrature T s and pressure p s for a surface luminosity L s are given by 

T s = (L s /4nR 2 a S B) 1/4 , (8) 

and 

p s = g s /K(p s ,T s ,Xi), (9) 

where g s = GM/R 2 , and <jsb — ac/4 is the Stefan-Boltzmann constant. The surface luminosity L s is 
determined so that the inner boundary condition L r = is satisfied at z = z . The envelope models 
calculated in this way are fully radiative with no convective regions, and are simply called radiative models 
in this paper. 

As discussed by Hanawa & Sugimoto (1982), there developes a convective zone in the nuclear burning 
region of the surface envelope immediately after the onset of a nuclear flash. To obtain a convective model 
in a nuclear flash, we first in the radiative models determine the location of the base of a convection zone 
as the inner most layer at which the condition for unstable helium burning 

(10) 

is satisfied (Cumming & Bildsten 2000), where e^ a denotes the 3a energy generation rate. We treat the 
mass extent z c — z e of the convection zone above the base as a parameter, where z c and z e are the mass 
depths at the base and at the top of the convection zone, respectively. The temperature gradient in the 
convection zone is given by equation (5), and the composition is determined by Xi = f* c dzXi(z)/ (z c — z e ), 
where Xi in the integrand is the composition in the corresponding radiative model. When we integrate the 
basic equations from the surface to the bottom of the envelope so that L r = at z = Zb, there appears a 
strong density discontinuity due to the discontinuous change of compositions Xj at the boundaries between 
the convection zone and the radiative regions, at which we assume the continuity of the pressure and the 
temperature. The envelope models thus constructed to have a convective zone are simply called convective 
models in this paper. 



de 3a 
dlnT 



^ r) max 
^ Ze CNO 



^ 1 { dhin 



4 V<91nT 



3. Oscillation Calculation 

For nonadiabatic and nonradial oscillations of rotating stars, we use the method of solution given 
by Lee & Saio (1987, 1993). Throughout this paper, we apply the Cowling approximation, neglecting 
the Euler perturbation of the gravitational potential, W. Since separation of variables is not possible 
for global oscillations of rotating stars, the perturbations are expanded in terms of spherical harmonic 
functions Y l m (9, <f>) for a given m, assuming axisymmetric equilibrium of the stars. The displacement vector 
£(r, 9, (j), t) is then given by 

oc 

tr = r Y, SiW^ m e i<rt , (11) 

2>|m| 

Z°=r E (^Yr + T v{ r)^YSr)e^ (12) 

l,V>\m\ V / 
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l,V>\m\ V Y 7 

and the Euler perturbation of the pressure, p'(r, 9, (f>, i), for example, is given by 

oo 

P'= E fi(0^ m e i<Tt , (14) 

;>|m| 

where <r is the oscillation frequency in an inertial frame, and I = \m\ + 2(j — 1) and I' = I + 1 for even 
modes, and Z = |m| + 2(j — 1) + 1 and I' = I — 1 for odd modes where j — 1, 2, 3, • • •. The expansions 
of pertubations are substituted into a set of linearized basic equations to give a set of simultaneous linear 
ordinary differential equations for the expansion coefficients Si(r), Hi(r), iTi'(r), p[(r), and so on. In a 
convective region, the total energy flux F is given as the sum of the radiative flux F ra d and the convective 
flux F conv so that F = F ra d + F conv . Since there exists no reliable theory for the interaction between 
pulsations and convective motion, for nonadiabatic calculations we assume for simplicity <5(V • F conv ) = 
and L ra( i = L r in the convective zone, where L r = 4irr 2 F r and L ra d — Aitr 2 F ra d y r- 

The inner mechanical boundary condtion we employ is given by p' = (McDermott & Taam 1987), 
and the inner thermal boundary condition is \Ss/c p \ « 1 (see, Lee & Baraffe 1995), where s is the specific 
entropy, and c p is the specific heat at constant pressure. The outer boundary conditions are Sp = and 
i5(47rr 2 o'sB7 l4 ) = SL r . The normalization condition max \iTi \ = 1 is imposed at the surface of the envelope. 
Because of the density discontinuity at the interfaces between the radiative regions and the convective zone, 
we need to impose jump conditions for the eigenfunctions, which are for adiabatic oscillations 

&.]t = 0, [Sp]t = 0, (15) 

where [/(x)]t = lim e _ +0 (/(a; + e) — f(x — e)). For nonadiabatic calculations, we need additional jump 
conditions, for which we employ for simplicity 

[SL rad }t=0, [5s]± = 0. (16) 

Truncation of the series expansions of perturbations is inevitable to get a set of linear ordinary 
differential equations of a finite dimesnion for numerical calculation. We find that including the first 
6 to 8 expansion coefficients is sufficient to obtain reasonbale convergence of the eigenfrequencies and 
eigenfunctions of the modes considered in this paper. 



4. Numerical Results 

4.1. Envelope Models 

For M = 1.4M Q , R = 10 6 cm, and z b /M Q = 1 x 1CT 10 and for the surface chemical composition X = 0.7, 
Z = 0.02, and Zcno = Z/2, we calculated envelope models for M = O.lM^dd and M = 0.02MEdd, where 
M E dd denotes the Eddignton accretion rate defined as M E dd = 4ncR/ K e = 1.88 x 10 18 (1 + X)" 1 (i?/10 6 )g 
s _1 and n e = 0.2(1 + X) is the electron scattering opacity. For the mass extent of the convection zone, we 
assume z e /z c = 0.01. In Table 1, we tabulate several physical quantities characterizing the envelopes, such 
as Ls/Lsdd, the radial thickness Ar of the envelope, the temperature Tb at the bottom of the envelope, 
and y c = z c /(At:R 2 ) for convective models, where L E dd = kncGM / K e is the Eddington luminosity. The 
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radial extent Ar and the base temperature of the convective envelope models are larger than those of 
the radiative models. The base z c of the convective zone resides in the mixed hydrogen-helium burning 
region for the model of M = OAMemi while the base is in the hydrogen-exhausted region for the model 
of M = 0.02Al Edd . In Figure 1, we plot the temperature T versus the column depth y = z/(AirR 2 ) for 
M = O.lMsdd in panel (a) and for M = 0.02MEdd in panel (b), where the solid lines and the dotted lines 
are for the radiative and the convective models, respectively. The temperature of the convective models is 
much higher than that of the radiative models. Both for the radiative and the convective models, the layers 
below the nuclear burning regions are nearly isothermal. To understand the modal properties of oscillation 
modes, it is useful to plot the Brunt- Vaisala frequency TV = y/—gA with g — GM r /r 2 and the Lamb 
frequency Li = + l)c s /r with c s being the sound velocity, where the Schwartzschild discriminant A is 
defined as 

_ d hip 1 rflnp pg 
dr Ti dr p 

where 

fd\np\ (dlnp\ fd\np\ d\nT dlnp 

and p is the mean molecular weight. In Figure 2, N 2 and Lf with I = 1 are plotted versus y for the 
models of M = O.lMEdd in panel (a) and for the models of M = 0.02MEdd in panel (b), where the solid 
lines and the dotted lines indicate respectively the radiative and the convective models, and iV 2 and Lf 
are normalized by GM/R 3 . Note that we have ./V 2 < in convective regions. Figure 3 gives the mean 
molecular weight p as a function of y, where the solid line and the dashed line are for the radiative and the 
convective models of M = OAMsdd, and the dotted line and the dash-dotted line are for the radiative and 
the convective models of M = 0.02MEdd- Figure 3 shows that the spiky or needle-like features found in 
the Brunt- Vaisala frequency is produced by extremely steep increase of the mean molecular weight p with 
increasing y in the active nuclear burning regions. 

4.2. Eigenfrequencies 

It is now well know that the corotating-frame oscillation frequency u> = a + mfl of the r-modes in 
non-isentropic stars is in the limit of ft — > asymptotic to the frequency given by 

w = 2mft/[J / (J / + 1)], (19) 

and that for a given combination of (m, V , f2) there exist the fundamental r-modc with no radial nodes of the 
eigenfunction %Ty and its overtones (see, e.g., Yoshida & Lee 2001). Note that the inertial-frame oscillation 



Table 1 . Physical parameters of the envelope models 



M/M Edd 


L s /L Edd 


Ar (cm) 


T b (K) 


Vc (g/cm 2 ) 


comment 


0.1 


2.70E-3 


3.17E3 


4.02E8 




radiative 




1.75E-2 


4.55E3 


9.70E8 


1.50E8 


convective 


0.02 


5.26E-4 


2.61E3 


2.01E8 




radiative 




4.34E-1 


6.82E3 


1.45E9 


1.65E8 


convective 



— (Vad- V)+V M 

Xp 



(17) 
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frcqucncy a of the r-modes is negative for positive m, indicating the modes are prograde observed in an 
inertial frame. In this paper, we calculate the fundamental r-modes of V = |m|, \m\ + 1, \m\ + 2, \m\ + 3 
for m — 1 and to = 2 at several values of fl for the envelope models. In Figure 4, the oscillation frequencies 
u)/2tt for the radiative models are plotted versus Q/2tt for M = 0.1M Edd in panel (a) and for M = 0.02M Edd 
in panel (b), where the solid lines with open symbols and those with filled symbols are for the modes of 
to = 1 and to — 2, respectively, and the circles, squares, triangles, and diamonds indicate the r-modes of 
/' = \m\, \m\ + 1, \m\ + 2, and \m\ + 3, respectively. As is increased, the frequency uu of the r-modes of 

V > \m\ deviates appreciably from the asymptotic frequency (19), and becomes almost insensitive to f2, 
gathering in a frequency range of lu/2tt <10Hz, at rapid rotation rates. For convenience, we have tabulated 
in Table 2 the oscillation frequency of the fundamental r-modes of to = 1 and 2 for 0/27T — 300 Hz and 
600 Hz for the radiaitve envelope models. Table 2 clearly shows that the r-mode frequency of the radiative 
models is also dependent on the mass accretion rates and is practically independent of the spin frequency at 
rapid rotation rates for I' > \m\. As shown in the next subsection, the amplitudes of the r-modes at rapid 
rotation rates are strongly confined in the equatorial regions of the star, in which the Coriolis parameter 
2fl sin 9 with being the latitude has small values. The modal properties found for the r-modes of a thin 
spherical shell at rapid rotation rates are quite similar to those of the low frequency equatorial waves in the 
ocean of the earth (e.g., Pedlosky 1987). 

Since in the propagation zone of the r-modes the pressure scale height H p « y/p of the convective 
models is larger than that of the radiative models because of higher envelope temperature of the former, 
the oscillation frequency u) of the r-modes in the convective models is increased compared to that in the 
radiative models (see the dispersion relation given in Discussion). This increase in to for the convective 
models in turn leads to a decrease in the inertial- frame oscillation frequency \a\, compared to that for 
the radiative models. In Figure 5, the relative difference / of the oscillation frequency of the fundamental 
r-modes between the radiative and the convective models is plotted versus a rad j2^ for M/Msdd = 0.1 in 
panel (a) and for M/Msdd = 0.02 in panel (b), where 

/ = -(fcomj - &rad)/(?rad 1 (20) 

and o CO nv and a rad are the oscillation frequencies for the convective and the radiative envelope models, 
respectively, and the solid lines with open symbols and those with filled symbols are for the r-modes of 
to = 1 and to = 2, respectively, and the circles, squares, triangles, and diamonds indicate the r-modes of 

V = |m|, |to| + 1, |to| + 2, and \m\ + 3, respectively. Note that the relative frequency difference / for the 
r-modes of (to, V) = (1, 1) is larger than 0.1, and is outside the range employed for the figure. As —<J rad is 
increased, / for the fundamental r-mode of given (to, V) decreases gradually after passing a maximum. We 
find that the relative differences / at \a vad \j2^ <~ 300Hz for m = 1 and at |oy a d|/27r ~ 600Hz for m = 2, for 
example, have similar values, and that / for given m and becomes smaller for larger values of I'. If we 
compare the relative differences / calculated here with the frequency drifts observed in burst oscillations 
in X-ray bursters, we find that the fundamental r-modes of I' > \m\ + 1 produce relative differences / 
consistent with observed frequency drifts at \<t\/2tt <~ 300Hz for to = 1 or at \<j\/2it <~ 600Hz for to = 2, 
although the differences / for the r-modes of I' — \m\ are too large to be reconciled with the observations 
for the models used in this paper. If the star is rotating as rapidly as Q/2ir <~ 600Hz, / for the to = 1 
r-modes also becomes consistent with observed values, which suggests the possibility that the burst QPO 
near 600Hz could be produced by the r-modes of to = 1 . We also note that to obtain values of / as small 
as ^0.001 wc need to consider the r-modes that have the index /' much larger than |m|. 

Although the r-modes in the radiative models are all pulsationally stable for the mass accretion rates 
considered in this paper (see Strohmayer & Lee 1996), the fundamental r modes in the convective models 



- 8 - 



are driven unstable by strong nuclear burning in the convective zone. In Figure 6, the growth timcscalc 
Tgrowth = — of the r-modes in the convective models is plotted as a function of Sl/2n for M = O.lMEdd 
in panel (a) and for M — 0.02MEdd m panel (b), where oj is the imaginary part of the eigenfrequency a, 
and the solid lines with open symbols and those with filled symbols are for the modes of m = 1 and m = 2, 
respectively, and the circles, squares, and triangles indicate the r-modes of /' = \m\, \m\ + 1, and \m\ + 2, 
respectively. T grow thS for M = 0.02Msdd are by about three order of magnitude shorter than those for 
M = O.lMEdd, because the envelope temperature of the former is much higher than that of the latter. Note 
that Tgrowth becomes almost insensitive to f2 for rapid rotation rates. 

4.3. Eigenfunctions 

The expansion coefficients iTy and Hi of the fundamental r-modcs of (to, /') = (2, 3) at 
n = n/^/GM/R 3 = 0.1 are shown versus the column depth y for the models with M = O.lMEdd 
in Figure 7, and those of (m, I') = (1, 1) for the models with M = 0.02Me<m in Figure 8, where the panels 
(a) and (b) are for the radiative models and the panels (c) and (d) for the convective models, and the 
amplitude normalization is given by max|iT;/| = 1 at the surface. As shown by Figure 7, the component 
iTi = \ m \ + i of the r-mode at Cl = 0.1 is no longer dominant over the other components of iTi> and Hi. This 
is also the case for the component i7j>_i m i of the r-mode of (to, I') = (1,1) shown in Figure 8. The modal 
properties of the modes at £1 = 0.1 arc thus quite different from the properties found for the r-modes in 
the limit of Ct —> 0, for which the toroidal component iTi> = \ m \ + i (or iTJ/ = | TO |) is dominant over the other 
components %Ty , Hi , and Si . This becomes clear if we show the 9 dependence of the eigenfunctions £e and 
£0, for example. To do so, we introduce the functions Xg(r, 9, <p) and X^r, 9, <p) defined as 

Xg(r, 6,4>)=^ [&(r, 9, <t>, t)e~ M /r] , (21) 

X^r, 9, (f>) - -3? [i^(r, 9, <j>, t)e~ i,7t /r] , (22) 

and, for r = R and <j> = 0, we plot the functions Xg (solid lines) and X^ (dotted lines) versus cos 9 for the 
fundamental r-modes of (to, I') = (1,1), (1,2), and (1,3) in Figure 9 and for those of (to, I') = (2,2), (2,3), 
and (2,4) in Figure 10, where the radiative model of M = 0.02Msdd has been used for the mode calculation, 
and the amplitude normalization is given by max|i7]'| = 1 at the surface. Here, in each figure the panels 
on the left-hand-side are for Vt = 0.01, and those on the right-hand-side are for f2 — 0.1. Although the 
r-modes have large amplitudes at the poles at f2 = 0.01, the amplitudes of the modes at Cl = 0.1 become 
well confined into the equatorial regions of the star, and have effectively no amplitudes at the poles. Note 
that the parity of the function Xgeg becomes the same as that of X^ with respect to cos 9. 

It may be instructive to give a plot of SL r for the r-modes considered in this paper, although 
nonadiabatic calculation is not necessarily quite reliable for the convective models because of our lack 
of knowledge concerning the interaction between pulsations and convection, and because of somewhat 
adhoc way of construction of the models. For the convective model of M / MEdd = 0.1, Figure 11 shows 
the expansion coefficients %i(6L r ^/L s ) as a function of y for the fundamental r-modc of (to,/') = (2,3) 
in panel (a) and for that of (to,/') = (2,5) in panel (b) at Q = 0.1, where the solid, dotted, dashed, and 
dash-dotted lines are for the expansion coefficients with / = \m\, \m\ + 2, \m\ + 4, and \m\ + 6, respectively, 
and the amplitude normalization is given by max|i7]/| = 1 at the surface. Figure 12 illustrates $t(8L r / L s ) 
versus cos#, where the solid and the dotted lines respectively designate the r-modes of (to,/') = (2,3) and 
(to, /') = (2, 5) at CI = 0.1, and we have assumed r = R and = 0. From Figures 11 and 12, we find that 
$i(5L r j/L s ) at the surface has large amplitudes, comparable to max|iT;/|, in the equatorial regions, and 
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that, although s R(SL r /L s ) of the r-modc of (m, I') = (2, 3) does not change its sign as a function of cos 6, 
$t(5L r /L s ) of the r-mode of (m,l') — (2,5) changes its sign at a middle latitude, which may result in a 
reduction of the detectability of the signals produced by the mode. Note that the r-modes of (m, I') = (2, 3) 
and (m,l') = (2,5) are even modes, and the perturbation SL r /L s is symmetric about the equator. On the 
other hand, the odd r-modes of I' = \m\, \m\ + 2, • • • have the eigenfunctions antisymmetric about the 
equator. 



5. Discussion 

Erecting at the equator a cartesian coordinate system in which the coordinates x, y, and z are directed 
to the east, the north, and the zenith, respectively, and applying the /3-plane approximation at the origin, 
we may obtain the dispersion relation for equatorial waves given by (Pedlosky 1987) 

Xlu 2 + * - P = (2j + 1)\/A, (23) 

where k is the dimensionless wavenumbcr in the x direction, and A is the separation constant between the 
horizontal and vertical momentum equations, and the velocity perturbation v' y , for example, have been 
assumed to have the form given by e l ( kx +^ 1 )^ (jj)G{z) . Here, the separation constant A is determined by 
the vertical-structure equation that governs the function G(z), and the non- negative-integer j comes from 
the normal mode equation for v' y , for which the function ^f(y) is shown to be proportional to e~ v ^ 2 Hj{ri), 
where Hj{rj) is the Hermite polynomial and 77 = V\y/L e (Pedlosky 1987). In the dispersion relation, the 
dimensionless oscillation frequency u) and the wavenumber k in the x direction are normalized as 

G)=w/(PoL e ), k = L e k, (24) 

where /?o = 20/i?, and 

L e = ^N D/f3 ( 25) 

is the equatorial Rossby internal-deformation radius, and A^ is the characteristic value of the Brunt- Vaisala 
frequency, D is the depth of the atmosphere or the ocean (Pedlosky 1987). Rewriting the dispersion relation 
using dimensional uj and k, and taking the low frequency limit, we obtain 

m k .Hfe = { 2 i + l)Vx, (20) 

which can be solved for u> to give 

U = ™ N0iD/R) , (27) 

(2j + 1)VA + N Dk 2 /f3 (2j + 1) VA + m 2 N (D/R)/2il 

where we have assumed k <~ m/R with m being an integer corresponding to the azumutal wave number. 
Note that the oscillation frequency is proportional to the depth D of the atmosphere or the ocean, and 
D/R << 1. In the limit of £1 — > 0, we obtain 



and we obtain at rapid rotation speeds 



u^2n/m, (28) 



^ruMD^ (2Q) 
(2j + l)VA 
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which does not depend on 17. If we regard D and N as the pressure-scale-height H p and the characteristic 
value of N found in the r-mode propagation region, the dispersion relation given above describes 
qualitatively well the fl dependence of the r-mode frequency u obtained in this paper. The quantity 
N (D/R) in equation (29) may depend on parameters such as the surface gravity g — GM/R 2 , the 
radius R, and the mass accretion rate, and extensive model calculations would be required to see how the 
quantity N (D/R) depends on the parameters. Note that the equatorial waves have large amplitudes at the 
equatorial regions but vanishing amplitudes at the poles as indicated by the function ^(y) oc e~ n ^ 2 Hj(rj), 
the propeties of which are similar to those of the r-modes found in this paper. 

6. Conclusion 

In this paper, we calculated the r-modes propagating in mass-accreting, nuclear burning, and 
geometrically thin envelopes on the surface of rotating neutron stars. In steady state approximation, we 
constructed thin envelope models which are fully radiative or have a convective region. The radiative and 
the convective models were intended to represent envelope structures in a late phase and in an early phase 
of an X-ray burst, respectively. The convective models have higher temperature and larger geometrical 
thickness than the corresponding radiative models. We find that the corotating- frame oscillation frequency 
uo of the fundamental r-modes of (to, I'), which is asymptotic tow= 2mCl/[l'(l' + 1)] when f2 — > 0, becomes 
essentially insensitive to f2 and gathers in a frequency range of uj/2tt <10Hz at rapid rotation rates. It is 
also found that the amplitudes of the r-modes at rapid rotation rates become well confined to the equatorial 
regions of the star. Although the fundamental r-modes of the radiative envelope models are all pulsationally 
stable for the mass accretion rates considered in this paper, the r-modes of the convective models are driven 
unstable by nuclear burning in the convection zone. 

It is found that the incrtial-frame oscillation frequency \a\ of the r-modes of given (m, Z',f2) for the 
convective models becomes lower than that for the corresponding radiative models. In order to estimate 
possible magnitudes of the relative frequency changes of the r-modes during an X-ray burst, we calculated 
/ = —{(Jconv — &rad) / &rad as a function of 0, where a r ad and a conv are the oscillation frequencies of the 
fundamental r-modes in the radiaitve and the convective envelope models, respectively. We find that, if 
we consider the burst oscillations of |cx ra ^ \/2ir <~ 300Hz or |<7 roe /|/27r ~ 600Hz, the fundamental r-modes 
of I' > |m| + 1 with to = 1 or to = 2 can produce the relative frequency changes / consistent with the 
observed relative frequency drifts of less than ~1%, although we need to consider the r-modes of /' much 
larger than |m| for values of / as small as ~0.001. We also find that / for the M — 0.02MEdd models 
is larger than that for the M = O.lAlsdd models, which may reflect the properties of nuclear flashes 
dependent on the mass accretion rates (see, e.g., Fujimoto, Hanawa, & Miyaji 1981) and could be tested 
observationally. Considering that the fundamental r-modes in the convective models are driven unstable 
by nuclear burning in the convective zone, it is tempting to postulate that the burst oscillations observed 
in X-ray bursts are produced by the r-modes of low indices (to, I') that are excited in the early phase of 
the nuclear flash and survive until the tail, changing the frequency slightly according to the changes in the 
envelope structure. The r-mode scenario is consistent with the observations that indicate the existence of 
two frequencies separated by a few Hz, since the fundamental r-modes of different Z"s for a given to could 
be excited simulataneously in a burst to have similar oscillation frequencies for rapidly rotating neutron 
stars. The temporal frequency decrease observed in the tail of bursts from 4U 1636-53 (Strohmayer 1999) 
and KS 1731-260 (Muno et al 2000) may be attributable to switching of the r-modes observable to us from 
those of {m,l'i) to (m,l' 2 ) where l[ > 1' 2 for frequency decrease. The mode switching will cause discontinuous 
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frequency changes and could be tested by a close observation of the frequency evolution. To make possible 
the close comparison between obervations and the r-mode scenario for burst oscillations, however, it is 
definitely necessary to use more realistic envelope models obtained from time-dependent flash calculation 
with a reliable treatment of the convection zone. 

Recently, Chakrabarty et al (2003) reported the detection of quasi-periodic oscillatons in X-ray bursts of 
the millisecond pulsar SAX J1808. 4-3658, showing that the oscillation frequencies agree quite well with the 
spin frequency of the star. The time evolution of the burst oscillation frequencies, however, was found quite 
different from that in other X-ray bursters that show burst oscillations. The frequency drift Av/v ~ 0.01 
is among the largest observed in any neutron star. The drift time scale is an order of magnitude faster 
than in the other neutron stars, and the oscillation frequency overshoots the spin frequency, reaching the 
maximum during the burst rise. The burst oscillations found in the pulsar are inconsistent with both the 
contracting shell model with angular momentum conservation and the r-mode scenario. A key factor for 
the burst oscillaitons in the star may be a strong magnetic field infered from the existence of coherent X-ray 
pulsations in non-burst phases, since the strong magnetic field has significant effects both on the modal 
properties of waves propagating in the surface layers and on the way of nuclear burning on the surface of 
the star. 
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Tabic 2. Oscillation frequency vf, = u> ra d/2Tr of the fundamental r- modes of the radiative envelopes 
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Fig. 1. — Temperature T as a function of the column depth y for M = O.lMEdd in panel (a) and for 
M = 0.02MEdd in panel (b), where the solid lines and the dotted lines arc for the radiative and the 
convective models, respectively. 
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Fig. 2. — Propagation diagram for the envelope models, where the square of Brunt- Vaisala, frequency N and 
Lamb frequency Li with I = 1 are displayed versus the column depth y for the radiative models (solid lines) 
and the convective models (dashed lines). TV 2 and Lf are normalized by GM/R 3 with M and R being the 
mass and the radius of the neutron star. Panels (a) and (b) are for the cases of M/MEdd = 0.1 and 0.02, 
respectively. 




Fig. 3. — Mean molecular weight \i as a function of the column depth y for the radiative (solid line) and 
the convective (dashsed line) models of M/MEdd — 0.1 and the radiative (dotted line) and the convective 
(dash-dotted line) models of M/M Edd = 0.02. 
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Fig. 4. — Oscillation frequencies u> of the fundamental r-modes of the radiative envelope models are given 
versus Q/2tt for M/MEdd = 0.1 in panel (a) and for M/MEdd = 0.02 in panel (b), where the solid lines 
attached by open symbols and those by filled symbols are for m = 1 and in = 2, respectively, and the circles, 
squares, triangles, and diamonds indicate the r-modcs of I' = \m\, \m\ + 1, \m\ + 2, \m\ + 3, respectively. 
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Fig. 5. — Relative frequency differences / = — {a conv — o ra d)lo ra d of the fundamental r-modes between the 
radiative and the convective envelope models are plotted versus arad/^ for M/MEdd = 0.1 in panel (a) 
and for M/MEdd = 0.02 in panel (b), where the solid lines attached by open symbols and those by filled 
symbols are for m = 1 and m = 2, respectively, and the circles, squares, triangles, and diamonds indicate 
the r-modes of /' = \m\, \m\ + 1, |m| + 2, |m| + 3, respectively. Here, a ra d and a conv are the inertial-frame 
oscillation frequencies of the r-modes in the radiaitve and the convective models, respectively. 
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Fig. 6. — Growth time-scales T growt h = — 1/07 of the fundamental r-modes of the convective envelope models 
are plotted versus fl/2ir for M/MEdd = 0.1 in panel (a) and for M/MEdd = 0.02 in panel (b), where the 
solid lines attached by open symbols and those by filled symbols are for m — 1 and m = 2, respectively, and 
the circles, squares, and triangles indicate the r-modes of = \m\, \m\ + 1, and \m\ + 2, respectively. Here, 
07 is the imaginary part of the cigcnfrcquency of the r-modes. 
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Fig. 7. — Expansion coefficients iTi> and Hi of the fundamental r-mode of (m,V) = (2,3) at O = 0.1 arc 
given versus the column depth y for M/MEdd = 0.1, where the solid, dotted, dashed, and dash-dotted lines 
indicate the coefficients associated with I = I' — 1 = \m\, \m\ + 2, \m\ + 4, and \m\ + 6, respectively, and 
the amplitude normalization is given by max|iXi'| = 1 at the surface. Here, panels (a) and (b) are for the 
radiative models, and panels (c) and (d) are for the convective models. 
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Fig. 8. — Expansion coefficients iTy and Hi of the fundamental r-mode of (to, Z') = (1, 1) at O = 0.1 are 
given versus the column depth y for M / MEdd — 0.02, where the solid, dotted, dashed, and dash-dotted lines 
indicate the coefficients associated with V — I — 1 = |to|, |to| + 2, \m\ + 4, and \m\ + 6, respectively, and 
the amplitude normalization is given by max|iT;/| = 1 at the surface. Here, panels (a) and (b) are for the 
radiative models, and panels (c) and (d) are for the convective models. 
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Fig. 9. — Functions Xe(r 7 9 7 (f>) and ^(r, for r = R and </> = are plotted versus cos# for the 
fundamental r-modcs of V = \m\, \m\ + 1, and \m\ + 2 for m = 1, where the solid lines and the dotted 
lines denote the functions Xg and X^, respectively, and the radiative model with M/MEdd — 0.02 has been 
used for the mode computation. The amplitude normalization is given by max \iTi> \ = 1 at the surface. 
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Fig. 11. — Expansion coefficients $t(5L r j/L s ) of the fundamental r-modes at Cl = 0.1 are given versus the 
column depth y for the convective model with M/Msdd = 0.1, where the solid, dotted, dashed, and dash- 
dotted lines are for the coefficients associated with I = |m|, \m\ + 2, \m\ +4, and \m\ + 6, respectively. Panel 
(a) is for the r-mode of (m,V) = (2,3), and panel (b) is for the r-mode of {m,V) = (2,5). The amplitude 
normalization is given by max |iT//| = 1 at the surface. 





Fig. 12. — ?R(5L r /L s ) at the surface is plotted as a function of cos# for the fundamental r-modes of 
(m,l') = (2,3) (solid line) and of (to, Z') = (2,5) (dotted line) at O = 0.1 for the convective model with 
M/MEdd = 0.1. The amplitude normalization is given by max|z7}/| = 1 at the surface. 



